*Code verified using Stata/SE 16.1.

//TABLE OF CONTENTS.
//0.0. Packages.
//1.0. Article.
//1.1. Some preparation.
//1.2. Descriptive Statistics.
//1.3. Regression Models.


//0.0. PACKAGES.

*We use "blindschemes" for all graphs, "covbal" for the balance statistics, and "xtabond2" for the dynamic GMM models. These packages can be downloaded via the Stata console:

ssc install blindschemes
ssc install covbal
ssc install xtabond2


//1.0. ARTICLE.

*Load dataset. Below code should be performed in order.


//1.1. SOME PREPARATION.

*Sort data and xtset
sort agency_id year
xtset agency_id year

*Generate inflation-adjusted budget.
gen i_budget=budget*inflation

*Generate log budget.
gen lni_budget=ln(i_budget)

*Generate budget growth.
bys agency_id: gen bgrowth = (D.i_budget/L.i_budget)
replace bgrowth=bgrowth*100

*Generate bootstrap identifier.
gen bs_agency_id=agency_id

*Use Hadi's ancient method to identify outliers.
hadimvo bgrowth, gen(odd)


//1.2. DESCRIPTIVE STATISTICS.

//Figure 1. Budget Growth by Year.

*Generate median growth per year.
by year, sort: egen mgrowth = median(bgrowth)

*Graph growth separated by government incumbent.
twoway (scatter bgrowth year if bgrowth>-100 & bgrowth<100 & incumbent==0, mcolor(gs12) msize(vsmall)) (scatter bgrowth year if bgrowth>-100 & bgrowth<100 & incumbent==1, mcolor(gs8) msymbol(triangle_hollow) msize(vsmall)) (line mgrowth year, lpattern(solid) lwidth(thin) lcolor(black)), scheme(plotplain) xtitle("Fiscal Year") ytitle("Budget Growth (%)", margin(0 0 0 0)) xtick(1970 (10) 2015) title("") legend(position(6) region(lcolor(none)) bmargin(zero) region(margin(zero)) rows(1) cols(2) order(1 "Social Democratic Gov't" 2 "Liberal-Conservative Gov't")) aspect(1) xsize(5) ysize(5) graphregion(m(l=0))


//Figure 2. Conflict Distribution by Year.

*Generate mean conflict per year.
by year, sort: egen conflict_tot=mean(conflict) if bgrowth~=.

*Graph conflict separated by government incumbent.
graph twoway (bar conflict_tot year if incumbent==0, lwidth(0) barw(0.9) fcolor(gs12)) (bar conflict_tot year if incumbent==1, lwidth(0) barw(0.9) fcolor(gs8)), scheme(plotplain) ylabel(0 (0.25) 1) xtitle(Fiscal Year) ytitle(Conflict Distribution, margin(0 0 0 0)) title("") ylabel(, format(%9.2fc)) xlabel(1970 (10) 2015) legend(position(6) region(lcolor(none)) bmargin(zero) region(margin(zero)) order(1 "Social Democratic Gov't" 2 "Liberal-Conservative Gov't") symxsize(2) rows(1) cols(2)) aspect(1) xsize(5) ysize(5) graphregion(m(l=0))


//1.3. REGRESSION MODELS.

*Sort data and xtset.
sort agency_id year
xtset agency_id year

*Note that we recode the outliers as missing values instead of restricting the estimation sample in the models. This is mainly to avoid some annoyances with the GMM functions in Model 3, which will otherwise incorporate information from the outlying values.
replace bgrowth=. if odd==1

*We also change the decimal format to avoid being burdened with unneccesary details.
set cformat %3.2f


//Table 1. Ideological Conflict and Agency Budget Growth in Sweden, 1971-2014.

*Model 1.
bootstrap, cluster(agency_id) idcluster(bs_agency_id) group(agency_id) seed(999) reps(1000): xtreg bgrowth i.conflict i.academic i.political i.private i.public i.council i.board i.creator cl.lni_budget i.year, fe

*Generate marginal predictions.
margins conflict

*Save marginal predictions for later use.
marginsplot, name(m1,replace) scheme(plotplain) plotopts(lcolor(black) lwidth(thin) mcolor(black) mlwidth(thin)) ciopts(lcolor(black) lwidth(medthin)) title("Model 1", size(medsmall) nospan) ylabel(0(1)5) yline(0 5, lpattern(dot) lwidth(thin) lcolor(gs10)) ytitle("Budget Growth (%)") xtitle("Ideological Conflict (0,1)") graphregion(m(l=0)) plotregion(m(l=0 r=0)) xscale(ra(-0.15 1.15))

*Model 2 (somewhat slow). 
bootstrap, cluster(agency_id) idcluster(bs_agency_id) group(agency_id) seed(999) reps(1000): mixed bgrowth i.conflict i.academic i.political i.private i.public i.council i.board i.creator cl.lni_budget i.year || agency_id: conflict, reml cov(un) iterate(100)
estat sd

*Save random-effects for later use. Note that we're saving the predicted marginal intercept, not the model intercept. 
predict re1*, reffects
gen slope_re1 = _b[1.conflict] + re11
margins conflict, post
gen intercept_re1 = _b[0bn.conflict] + re12
gen yhat_re1= intercept_re1 + (slope_re1*conflict)

*Save marginal predictions for later use.
marginsplot, name(m2,replace) scheme(plotplain) plotopts(lcolor(black) lwidth(thin) mcolor(black) mlwidth(thin)) ciopts(lcolor(black) lwidth(medthin)) title("Model 2", size(medsmall) nospan) ylabel(0(1)5) yline(0 5, lpattern(dot) lwidth(thin) lcolor(gs10)) ytitle("Budget Growth (%)") xtitle("Ideological Conflict (0,1)") graphregion(m(l=0)) plotregion(m(l=0 r=0)) xscale(ra(-0.15 1.15))

*Model 3.
xtabond2 bgrowth l.bgrowth i.conflict academic political private public council board creator l.lni_budget y1974-y2014, gmmstyle(l.bgrowth conflict academic political private public council board l.lni_budget, collapse l(1 4)) ivstyle(creator y1974-y2014, eq(level)) robust twostep orthogonal

*Generate marginal predictions.
margins conflict

*Save marginal predictions for later use.
marginsplot, name(m3,replace) scheme(plotplain) plotopts(lcolor(black) lwidth(thin) mcolor(black) mlwidth(thin)) ciopts(lcolor(black) lwidth(medthin)) title("Model 3", size(medsmall) nospan) ylabel(0(1)5) yline(0 5, lpattern(dot) lwidth(thin) lcolor(gs10)) ytitle("Budget Growth (%)") xtitle("Ideological Conflict (0,1)") graphregion(m(l=0)) plotregion(m(l=0 r=0)) xscale(ra(-0.15 1.15))


//Figure 3. Predicted Marginal Budget Growth.

*Combine marginal predictions from Model 1, 2, and 3.
gr combine m1 m2 m3, ycommon title("") scheme(plotplain) row(1) imargin(0 0 0 0) ysize(2) scale(2)


//Figure 4. Predicted Agency-specific Budget Growth.

*To clean up the agency-specific graph, we only include agencies with actual variance on the treatment variable. This requires constructing measures that identify whether a given agency features both 0s and 1s on conflict.
sort agency_id conflict
egen xmin=min(conflict), by(agency_id)
egen xmax=max(conflict), by(agency_id)
gsort -xmin -xmax agency_id conflict

*We also want to tag the marginal predictions for reference.
gen margins_re1=3.70 if conflict==0
replace margins_re1=1.74 if conflict==1

*Graph the random-effects we saved above.
twoway (connected yhat_re1 conflict if xmin==0 & xmax==1, connect(L) lwidth(thin) lpattern(solid) lcolor(gs12) mstyle(none)) (connected margins_re1 conflict, sort lwidth(thin) lpattern(solid) lcolor(black) msymbol(circle_hollow) mlcolor(black) mlwidth(thin)) (connected yhat_re1 conflict if xmin==0 & xmax==1 & agency_id==9089, sort lwidth(thin) lpattern(longdash) lcolor(black) msymbol(circle_hollow) mlcolor(black) mlwidth(thin)) (connected yhat_re1 conflict if xmin==0 & xmax==1 & agency_id==5002, sort lwidth(thin) lpattern(shortdash) lcolor(black) msymbol(circle_hollow) mlcolor(black) mlwidth(thin)), name(re1,replace) title("") ylabel(10(2)-6) yscale(range(-6.5 10.5)) xscale(range(-0.15 1.15)) ytitle("Budget Growth (%)") xlabel(0 1) xtitle("Ideological Conflict (0,1)") scheme(plotplain) legend(off) r2title(" " " " " ")  text(1.74 1 "  Average Marginal" " Effect", size(small) place(e)) text(7.6 1 "  Ethnic Discrimination" " Ombudsman", size(small) place(e)) text(-4.8 1 "  Swedish State" "  Railways", size(small) place(e)) aspect(1) xsize(6) ysize(5) graphregion(m(l=0 r=0)) plotr(m(zero))



